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ABSTRACT 

In  this  chapter,  we  consider  the  sensor  management  problem  for  tracking 
space  targets  where  the  targets  may  apply  evasive  maneuvering  strategy  to  avoid 
being  tracked  by  the  space  borne  observers.  We  first  study  the  case  of  single  target 
tracking  by  a  single  observer  and  formulate  the  pursuit-evasion  game  with  complete 
information.  Then  we  extend  the  tracking  problem  to  a  set  of  collaborative  observers 
and  each  observer  has  to  decide  when  to  sense  which  target  in  order  to  achieve  the 
desired  estimation  error  covariance.  A  popularly  used  criterion  for  sensor  management 
is  to  maximize  the  total  information  gain  in  the  observer-to-target  assignment. 

We  compare  the  information  based  approach  to  the  game  theoretic  criterion  where 
the  observers  are  assigned  according  to  the  best  response  of  the  terminal  result  in 
the  pursuit— evasion  game.  Finally,  we  use  realistic  satellite  orbits  to  simulate  the 
space  resource  management  for  situation  awareness.  We  adopted  NASA’s  General 
Mission  Analysis  Tool  (GMAT)  for  space  target  tracking  with  multiple  space  borne 
observers.  The  results  indicate  that  the  game  theoretic  approach  is  more  effective 
than  the  information  based  approach  in  handling  intelligent  target  maneuvers. 


Chapter  12 

Orbital  Evasive  Target  Tracking  and  Sensor 
Management 


Huimin  Chen,  Genshe  Chen,  Dan  Shen, 
Erik  P.  Blasch,  and  Khanh  Pham 


Summary  In  this  chapter,  we  consider  the  sensor  management  problem  for  track¬ 
ing  space  targets  where  the  targets  may  apply  evasive  maneuvering  strategy  to  avoid 
being  tracked  by  the  space  borne  observers.  We  first  study  the  case  of  single  target 
tracking  by  a  single  observer  and  formulate  the  pursuit-evasion  game  with  complete 
information.  Then  we  extend  the  tracking  problem  to  a  set  of  collaborative  observers 
and  each  observer  has  to  decide  when  to  sense  which  target  in  order  to  achieve  the 
desired  estimation  error  covariance.  A  popularly  used  criterion  for  sensor  manage¬ 
ment  is  to  maximize  the  total  information  gain  in  the  observer-to-target  assignment. 
We  compare  the  information  based  approach  to  the  game  theoretic  criterion  where 
the  observers  are  assigned  according  to  the  best  response  of  the  terminal  result  in 
the  pursuit-evasion  game.  Finally,  we  use  realistic  satellite  orbits  to  simulate  the 
space  resource  management  for  situation  awareness.  We  adopted  NASA’s  General 
Mission  Analysis  Tool  (GMAT)  for  space  target  tracking  with  multiple  space  borne 
observers.  The  results  indicate  that  the  game  theoretic  approach  is  more  effective 
than  the  information  based  approach  in  handling  intelligent  target  maneuvers. 


12.1  Introduction 

Over  recent  decades,  the  space  environment  has  become  more  complex  with  a  sig¬ 
nificant  increase  in  space  debris  among  densely  populated  satellites.  Efficient  and 
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reliable  space  operations  rely  heavily  on  the  space  situation  awareness  where  search 
and  tracking  space  targets  and  identifying  their  intent  are  crucial  in  creating  a  con¬ 
sistent  global  picture  of  the  space.  Orbit  determination  with  measurements  pro¬ 
vided  by  a  constellation  of  satellites  has  been  studied  extensively  [6,  7,  14,  16]. 
The  tracking  and  data  relay  satellite  system  uses  satellites  in  geostationary  orbits 
(GEO)  to  track  the  satellites  in  low-Earth  orbit  (LEO).  The  global  positioning  sys¬ 
tem  (GPS)  uses  a  constellation  of  satellites  with  pseudo-range  measurements,  i.e., 
range  measurements  with  clock  differentials  to  determine  the  location  of  a  user. 
Unlike  ground  targets  whose  motion  may  contain  frequent  maneuvers,  a  satellite 
usually  follows  its  orbit  so  that  long-term  prediction  of  its  orbital  trajectory  is  pos¬ 
sible  once  the  orbital  elements  are  known  [6].  However,  a  space  target  can  also 
make  an  orbital  change  owing  to  its  desired  mission  or  intentionally  hiding  from  the 
space  borne  observers.  Existing  maneuvering  target  tracking  literature  mainly  fo¬ 
cuses  on  modeling  target  maneuver  motion  at  random  onset  time  (see,  e.g.,  [1,  13]). 
In  space  surveillance,  the  constellation  of  satellite  observers  is  usually  known  to 
the  adversary  and  very  unlikely  to  change  frequently  due  to  energy  constraint.  In 
this  case,  evasive  maneuvering  motion  can  be  intelligently  designed  to  take  advan¬ 
tage  of  the  sensing  geometry,  e.g.,  transferring  to  an  orbit  with  maximum  dura¬ 
tion  of  the  Earth  blockage  to  an  observer  with  known  orbital  trajectory.  Accord¬ 
ingly,  a  sensor  management  method  has  to  optimally  utilize  the  sensing  resources 
to  acquire  and  track  space  targets  with  sparse  measurements,  i.e.,  with  a  typically 
large  sampling  interval,  in  order  to  maintain  a  large  number  of  tracks  simultane¬ 
ously. 

Sensor  management  is  concerned  with  the  sensor-to-target  assignment  and  a 
schedule  of  sensing  actions  for  each  sensor  in  the  near  future  given  the  currently 
available  information  on  the  space  targets.  Sensor  assignment  and  scheduling  usu¬ 
ally  aim  to  optimize  a  certain  criterion  under  energy,  data  processing  and  communi¬ 
cation  constraints.  One  popularly  used  criterion  is  the  total  information  gain  for  all 
the  targets  being  tracked  [11].  However,  this  criterion  does  not  prioritize  the  targets 
with  respect  to  their  types  or  identities.  Alternatively,  covariance  control  optimizes 
the  sensing  resources  to  achieve  the  desired  estimation  error  covariance  for  each 
target  [10].  It  has  the  flexibility  to  design  the  desired  tracking  accuracy  according  to 
the  importance  of  each  target.  We  want  to  compare  both  informative  based  criterion 
and  the  covariance  control  method  in  a  game  theoretic  setting  where  the  target  can 
intelligently  choose  its  maneuvering  motion  and  onset  time  based  on  the  knowledge 
of  observers’  constellation. 

The  rest  of  the  chapter  is  organized  as  follows.  Section  12.2  presents  the 
space  target  motion  and  sensor  measurement  model.  Section  12.3  provides  a  game 
theoretic  formulation  of  target  maneuvering  motion.  Section  12.4  discusses  the 
information  based  performance  metric  for  sensor  management  and  covariance 
control  method  with  possibly  evasive  target  maneuvering  motion.  Section  12.5 
compares  the  performance  of  the  proposed  tracking  and  sensor  management 
scheme  with  the  existing  methods.  Conclusions  and  future  work  are  presented  in 
Sect.  12.6. 
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12.2  Fundamentals  of  Space  Target  Orbits 
12.2.1  Time  and  Coordinate  Systems 

Several  time  systems  are  popularly  used  in  the  orbit  determination  problems.  Satel¬ 
lite  laser  ranging  measurements  are  usually  time-tagged  in  coordinated  universal 
time  (UTC)  while  global  positioning  system  (GPS)  measurements  are  time  tagged 
in  GPS  system  time  (GPS-ST).  Although  both  UTC  and  GPS-ST  are  based  on 
atomic  time  standards,  UTC  is  loosely  tied  to  the  rotation  of  the  Earth  through  the 
application  of  “leap  seconds”  while  GPS-ST  is  continuous  with  the  relation  GPS- 
ST  =  UTC  +  n  where  n  is  the  number  of  leap  seconds  since  January  6,  1980.  The 
orbital  equation  describing  near-Earth  satellite  motion  is  typically  tagged  with  ter¬ 
restrial  dynamical  time  (TDT).  It  is  an  abstract,  uniform  time  scale  implicitly  de¬ 
fined  by  the  motion  equation  and  can  be  converted  to  UTC  or  GPS-ST  for  any  given 
reference  date. 

The  Earth  centered  inertial  (ECI)  coordinate  system  used  to  link  GPS-ST  with 
UTC  is  a  geocentric  system  defined  by  the  mean  equator  and  vernal  equinox  at 
Julian  epoch  2000.0.  Its  XT-plane  coincides  with  the  equatorial  plane  of  the  Earth 
and  the  X-axis  points  toward  the  vernal  equinox  direction.  The  Z-axis  points  toward 
the  north  pole  and  the  T-axis  completes  the  right  hand  coordinate  systems. 

The  Earth  centered  Earth  fixed  (ECEF)  coordinate  system  has  the  same  X  T  -plane 
and  the  Z-axis  as  in  the  inertial  coordinate  system.  However,  its  X-axis  rotates  with 
the  Earth  and  points  to  the  prime  meridian  and  the  Y -axis  completes  the  right  hand 
coordinate  systems. 

The  local  Cartesian  system  commonly  referred  to  as  east-north-up  (ENU)  coor¬ 
dinate  system  has  its  origin  at  some  point  on  the  Earth  surface  or  above  (typically  at 
the  location  of  an  observer).  Its  Z-axis  is  normal  to  the  Earth’s  reference  ellipsoid 
defined  by  the  geodetic  latitude.  The  X-axis  points  toward  the  east  while  the  T-axis 
points  toward  the  north.  The  conversion  among  these  three  coordinate  systems  is 
provided  in  Appendix  1 . 


12.2.2  Orbital  Equation  and  Orbital  Parameter  Estimation 


Without  any  perturbing  force,  the  position  r  of  a  space  target  relative  to  the  center 
of  the  Earth  in  ECI  coordinate  system  should  satisfy 


r  = 


(12.1) 


where  //  is  the  Earth’s  gravitational  parameter.  The  target  velocity  is  v  =  r  and  the 
radial  velocity  is  vr  =  ™  where  r  =  ||r||  is  the  distance  from  the  target  to  the  cen¬ 
ter  of  the  Earth.  In  order  to  determine  the  position  and  velocity  of  a  satellite  at  any 
time  instance,  six  parameters  are  needed,  typically,  the  three  position  components 
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and  three  velocity  components  at  a  certain  time  instance.  Alternatively,  the  orbital 
trajectory  can  be  conveniently  described  by  the  six  components  of  the  Keplerian  ele¬ 
ments.  The  description  of  Keplerian  elements  and  their  relationship  to  the  kinematic 
state  of  the  target  can  be  found  in  Appendix  2. 

In  reality,  a  number  of  forces  act  on  the  satellite  in  addition  to  the  Earth’s  gravity. 
To  distinguish  them  from  the  central  force  created  by  the  satellite  target,  these  forces 
are  often  referred  to  as  perturbing  forces.  In  a  continuous  time  state  space  model, 
perturbing  forces  are  often  lumped  into  the  noise  term  of  the  system  dynamics. 
Denote  by  x(f)  the  continuous  time  target  state  given  by 


x(f)  = 


r(f) 

m 


x(t) 

y(t ) 

z(t) 

1 hit) 

Vy(t) 

VZ(0 


(12.2) 


For  convenience,  we  omit  the  argument  t  and  write  the  nonlinear  state  equation  as 
follows. 


x=  /(x)  +  W 


(12.3) 


where 


/(x)  = 


Vx 

Vy 

VZ 

-Or/r3)x 

~(n/r3)y 

~(m/ f3)z 


(12.4) 


and 


0 

0 

0 

wx 

Wy 

wz 


(12.5) 


is  the  acceleration  resulting  from  perturbing  forces.  As  opposed  to  treating  the  per¬ 
turbing  acceleration  as  noise,  spacecraft  general  propagation  (SGP)  model  maintains 
general  perturbation  element  sets  and  finds  analytical  solution  to  the  satellite  motion 
equation  with  time  varying  Keplerian  elements  [12],  For  precise  orbit  determination, 
numerical  integration  of  (12.3)  is  often  a  viable  solution  where  both  the  epoch  state 
and  the  force  model  have  to  be  periodically  updated  when  a  new  measurement  is 
available  [17]. 
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12.3  Modeling  Maneuvering  Target  Motion  in  Space  Target 
Tracking 

12.3.1  Sensor  Measurement  Model 

We  consider  the  case  that  a  space  satellite  in  low-Earth  orbit  (LEO)  observes  a 
target  in  geostationary  orbit  (GEO).  A  radar  onboard  the  space  satellite  can  provide 
the  following  type  of  measurements:  range,  azimuth,  elevation,  and  range  rate.  The 
range  between  the  ith  observer  located  at  (x,- ,  y,-,  Zi)  and  the  space  target  located  at 
( x,y,z )  is  given  by 


dr(i )  =  J(x-  X()2  +  (y  -  v/)2  +  (s  -  Zi)2 


(12.6) 


The  azimuth  is 


The  elevation  is 


da(i )  =  tan 


-l  (y-  yt 

X  -  Xj 


de(i)  =  tan  1  ( 

V-y/(x  —  Xj)2  +  (y  —  Vi)2 

The  range  rate  is 

(x  -  x,)(x  -  xj)  +  (y  -  y;)(y  -  y,)  +  (z  -  z,)(z  -  Zi) 
drU)  = - ; - 


(12.7) 


(12.8) 


(12.9) 


Measurements  from  the  th  observer  will  be  unavailable  when  the  line-of-sight 
path  between  the  observer  and  the  target  is  blocked  by  the  Earth.  Thus,  the  con¬ 
stellation  of  multiple  observers  is  important  to  maintain  consistent  coverage  of  the 
target  of  interest. 

The  condition  of  Earth  blockage  is  examined  as  follows.  If  there  exist  a  e  [0,  1] 
such  that  Da{i)  <  Re,  where 


Da{i)  =  \/[(l  -  «)x,-  +  ax]2  +  [(1  -  ar)y;  +  ay]"  +  [(1  -  a)zi  +  az\~  (12. 


10) 


then  the  measurement  from  the  i  th  observer  to  the  target  will  be  unavailable.  The 
minimum  of  Da(i)  is  achieved  at  a  =  a*  given  by 


x,  (x  -  Xj)  +  y,-(y  -  y,-)  +  Zj(z  -  Zj) 
(x  -  x,-)2  +  (y  -  yi)2  +  (z-  Zi)2 


(12.11) 


Thus,  we  first  examine  whether  a*  e  [0,  1]  and  then  check  the  Earth  blockage  con¬ 
dition  Da*(i)  <  Re- 


238 


H.  Chen  et  al. 


12.3.2  Game  Theoretic  Formulation  for  Target  Maneuvering 
Onset  Time 

We  consider  the  case  that  a  single  observer  tracks  a  single  space  target.  Initially, 
the  observer  knows  the  target’s  state  and  the  target  also  knows  the  observer’s  state. 
Assume  that  the  target  can  only  apply  a  T  second  burn  that  produces  a  specific  thrust 
w  with  a  maximum  acceleration  of  a  m/s2.  The  goal  of  the  target  is  to  determine 
the  maneuvering  onset  time  and  the  direction  of  the  thrust  so  that  the  resulting  orbit 
will  have  the  maximum  duration  of  the  Earth  blockage  to  the  observer.  The  goal  of 
the  observer  is  to  maintain  the  target  track  with  the  highest  estimation  accuracy.  To 
achieve  this,  the  observer  has  to  determine  the  sensor  revisit  time  and  cuing  region  as 
well  as  notify  other  observers  having  better  geometry  when  Earth  blockage  occurs. 
Without  loss  of  generality,  we  assume  that  the  target  can  transfer  its  orbit  to  the 
same  plane  as  the  observer.  In  this  case,  when  the  target  is  at  the  opposite  side  of  the 
Earth  with  respect  to  the  observer  and  rotating  in  the  same  direction  as  the  observer, 
the  duration  of  the  Earth  blockage  will  be  the  maximum  compared  with  other  orbits 
with  the  same  orbital  elements  except  the  inclinations.  Note  that  in  the  pursuit- 
evasion  game  confined  to  a  two  dimensional  plane,  the  minimax  solution  requires 
that  the  target  applies  the  same  thrust  angle  as  the  observer’s  (see  Appendix  4  for 
details).  Thus,  an  intelligent  target  will  choose  its  maneuvering  onset  time  as  soon  as 
its  predicted  observer’s  orbital  trajectory  has  the  Earth  blockage.  The  corresponding 
maneuvering  thrust  will  follow  the  minimax  solution  to  the  pursuit-evasion  game. 

When  a  target  is  tracked  by  multiple  space  borne  observers,  an  observer  can 
predict  the  target’s  maneuvering  motion  based  on  its  estimated  target  state  and  the 
corresponding  response  of  the  pursuit-evasion  game  from  the  target  where  the  ter¬ 
minal  condition  will  lead  to  the  Earth  blockage  to  the  observer.  Thus,  there  is  a  need 
for  the  sensor  manager  to  select  the  appropriate  set  of  sensors  that  can  persistently 
monitor  all  the  targets  especially  when  they  maneuver. 


12.3.3  Nonlinear  Filter  Design  for  Space  Target  Tracking 

When  a  space  target  has  been  detected,  the  filter  will  predict  the  target  state  at  any 
time  instance  in  the  future  based  on  the  available  sensor  measurements.  Denote  by 
xjT  the  state  prediction  from  time  t\-  \  to  time  4  based  on  the  state  estimate  xjf_  j  at 
time  4_ 1  with  all  measurements  up  to  4_i .  The  prediction  is  made  by  numerically 
integrating  the  state  equation  given  by 

m  =  f(ut))  (i2.i2) 

without  process  noise.  The  mean  square  error  (MSE)  of  the  state  prediction  is  ob¬ 
tained  by  numerically  integrating  the  following  matrix  equation 

Ht)  =  F{x^)P(t)  +  P{t)F(x~)T  +  Q{t ) 


(12.13) 
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where  F(x 7)  is  the  Jacobian  matrix  given  by 

F(x)  = 


03x3  h 
Fo{x)  03x3 
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(12.14) 


(12.15) 


(12.16) 


and  evaluated  at  x  =  xk  .  The  measurement  zk  obtained  at  time  tk  is  given  by 


zk  —  h(xk)  +  \k 


(12.17) 


where 

vk~AT(0,Rk)  (12.18) 

is  the  measurement  noise,  which  is  assumed  independent  of  each  other  and  indepen¬ 
dent  to  the  initial  state  as  well  as  process  noise. 

The  recursive  linear  minimum  mean  square  error  (LMMSE)  filter  applies  the 
following  update  equation  [2] 


*k\k  =  E*[xk  |ZA  ]  =  Xk\k-l  +  KkZk\k-l  (12.19) 

Pk\k  =  Pk\k-1  -  KkSkK'k  (12.20) 


where 


xk\k-i  =  E*[xk\Zk~l] 
zk\k-i  =  E*[zk\Zk~l] 

*-k\k- 1  =  X*  —  XjtlJ;-! 

h\ k- 1  =  Ik  -  %k\k- 1 

pk\k- 1  =  E[h\k-\xk\k-\\ 

Sk  =  E[zk\k-iz'k^_i] 

Ek  =  Ckkik  Sk 
C\kzk  —  E\xk \k-\z!k\k_\\ 

Note  that  £*[•]  becomes  the  conditional  mean  of  the  state  for  linear  Gaussian  dy¬ 
namics  and  the  above  filtering  equations  become  the  celebrated  Kalman  filter  [2]. 
For  nonlinear  dynamic  system,  (12.19)  is  optimal  in  the  mean  square  error  sense 
when  the  state  estimate  is  constrained  to  be  an  affine  function  of  the  measurement. 
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Given  the  state  estimate  x.k-\\k-t  and  its  error  covariance  1|*_  1  at  time  tk—\, 

if  the  state  prediction  Xk\k-i,  the  corresponding  error  covariance  Pk\k- 1,  the  mea¬ 
surement  prediction  %k\k-i ,  the  corresponding  error  covariance  Sk,  and  the  crossco¬ 
variance  £[xj:|/t-iz[.|/(,_1]  in  (12.19)  and  (12.20)  can  be  expressed  as  a  function  only 
through  x*—  1 1*— i  and  Pk—\\k-i,  then  the  above  formula  is  truly  recursive.  However, 
for  general  nonlinear  system  dynamics  (12.3)  and  measurement  equation  (12.17)  , 
we  have 


X-k\k-l 


/(x(f),  w(t))dt  +  xk-l\Zk  1 


z  k\k-i  =  E*[h(jLk,Vk)\Zk-1] 


(12.21) 

(12.22) 


Both  Xk\k- 1  and  ik\k- 1  will  depend  on  the  measurement  history  Zk~ 1  and  the  corre¬ 
sponding  moments  in  the  LMMSE  formula.  In  order  to  have  a  truly  recursive  filter, 
the  required  terms  at  time  4  can  be  obtained  approximately  through  x* _ i|*— t  and 

Pk-I\k-1,  i-e., 


{x*|*-l>  Pk\k-l}  ~  Pred[/ (■),  x*_i|*_i ,  Pk-\\k-\] 

{^k\k—  1;  ^  Pred  [/?(■),  Xk\k-l,  Pk\k— l] 

where  Pred[/(-),  x*-i|*-i,  Pk-i\k-l]  denotes  that  {x*_i|*_!,  Pk_X\k-\}  propagates 
through  the  nonlinear  function  /(•)  to  approximate  £*[/(•) |Zi_l]  and  the  corre¬ 
sponding  error  covariance  Pk\k-\  ■ 

Similarly,  Pred[/?(-),  Xk\k-i ,  Pk\k-i]  predicts  the  measurement  and  the  corre¬ 
sponding  error  covariance  only  through  the  approximated  state  prediction.  This 
poses  difficulties  for  the  implementation  of  the  recursive  LMMSE  filter  due  to  insuf¬ 
ficient  information.  The  prediction  of  a  random  variable  going  through  a  nonlinear 
function,  most  often,  can  not  be  completely  determined  using  only  the  first  and  sec¬ 
ond  moments.  Two  remedies  are  often  used:  One  is  to  approximate  the  system  to  the 
best  extent  such  that  the  prediction  based  on  the  approximated  system  can  be  carried 
out  only  through  {x^-ii/t-i,  Pk-\\k-\)  [20].  Another  is  by  approximating  the  den¬ 
sity  function  with  a  set  of  particles  and  propagating  those  particles  in  the  recursive 
Bayesian  filtering  framework,  i.e.,  using  a  particle  filter  [8], 


12.3.4  Posterior  Cramer-Rao  Lower  Bound  of  the  State 
Estimation  Error 

Denote  by  J  (t)  the  Fisher  information  matrix.  Then  the  posterior  Cramer-Rao  lower 
bound  (PCRLB)  is  given  by  [18] 


-l 


B(t)  =  J(t ) 


(12.23) 
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which  quantifies  the  ideal  mean  square  error  of  any  filtering  algorithm,  i.e., 

£[(x(4)  -  x(f;fc)) (x(4)  -  x(tk))T\Zk]  >  B(tk)  (12.24) 

Assuming  an  additive  white  Gaussian  process  noise  model,  the  Fisher  information 
matrix  satisfies  the  following  differential  equation 

Ht)  =  (x)  -  F(x)TJ{t)  -  J (t)Q(t)J (?)  (12.25) 


for  tk- 1  <  t  <  tk  where  F  is  the  Jacobian  matrix  given  by 

3/(x) 


Fix)  = 


dx 


(12.26) 


When  a  measurement  is  obtained  at  time  tk  with  additive  Gaussian  noise  A/"(0,  Rk), 
the  new  Fisher  information  matrix  is 


J(t+)  =  J(t~)  +  Ex[H(x)TR~lH(x )] 
where  H  is  the  Jacobian  matrix  given  by 


(12.27) 


H(x)  = 


3ft  (x) 
dx 


(12.28) 


See  Appendix  3  for  the  numerical  procedure  to  evaluate  the  Jacobian  matrix  for  a 
non-perturbed  orbital  trajectory  propagation.  The  initial  condition  for  the  recursion 
is  J(t0)  and  the  PCRLB  can  be  obtained  with  respect  to  the  true  distribution  of 
the  state  x(f).  In  practice,  the  sensor  manager  will  use  the  estimated  target  state  to 
compute  the  PCRLB  for  any  time  instance  of  interest  and  decide  whether  a  new 
measurement  has  to  be  made  to  improve  the  estimation  accuracy. 


12.4  Sensor  Management  for  Situation  Awareness 

12.4.1  Information  Theoretic  Measure  for  Sensor  Assignment 

In  sensor  management,  each  observer  has  to  decide  when  to  measure  which  target 
so  that  the  performance  gain  in  terms  of  a  certain  metric  can  be  maximized.  For 
a  Kalman  filter  or  its  extension  for  the  nonlinear  dynamic  state  or  measurement 
equations,  namely,  the  recursive  LMMSE  filter,  the  error  covariance  of  the  state 
estimate  has  the  following  recursive  form  [2]. 

n+V+t  =  Pf+x\k  +  H{xk+l)T  RklxH{xk+l)  (12.29) 

Thus  the  information  gain  from  the  sensor  measurement  at  time  tk+\  in  terms  of  the 
inverse  of  the  state  estimation  error  covariance  is  H(xk+ 1)  Rk+[H (xk+\)  where 
Rk+ 1  is  the  measurement  error  covariance.  Consider  M  observers  each  of  which  can 
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measure  at  most  one  target  at  any  sampling  time.  When  there  are  N  space  targets 
being  tracked  by  M  observers,  sensor  assignment  is  concerned  with  the  sensor  to 
target  correspondence  so  that  the  total  information  gain  can  be  maximized.  Denote 
by  Xij  the  assignment  of  observer  i  to  target  j  at  any  particular  time  4+ 1  •  The  sensor 
assignment  problem  is 


mm  a  j  Xij 
Xij 

M 

(12.30) 

subject  to 

J2  Kj  -  1  ’  2  =  1,- 
!  =  1 

N 

,.,N- 

(12.31) 

i=l>- 
j= 1 

.  ,,M; 

(12.32) 

and  Xij  £  { 0 ,  1}.  The  cost  Cjj  is 

Cij=Ti{H (xj ((4+ 1 )) ) 7  R~ 1  (4+t ) ^  (xj  ( (4+i)) ) }  (12.33) 

if  there  is  no  Earth  blockage  between  observer  i  to  target  j .  In  the  above  formulation, 
all  targets  are  assumed  to  have  the  same  importance  so  the  observes  are  scheduled  to 
make  the  most  informative  measurements.  This  may  lead  to  a  greedy  solution  where 
those  targets  close  to  the  observers  will  be  tracked  more  accurately  than  those  away 
from  the  observers.  It  may  not  achieve  the  desired  tracking  accuracy  for  each  target. 


12.4.2  Covariance  Control  for  Sensor  Scheduling 

Covariance  control  method  does  not  optimize  a  performance  metric  directly  in  the 
sensor-to-target  assignment,  instead,  it  requires  the  filter  designer  to  specify  a  de¬ 
sired  state  estimation  error  covariance  so  that  the  selection  of  sensors  and  the  corre¬ 
sponding  sensing  times  will  meet  the  specified  requirements  after  the  tracker  update 
using  the  sensor  measurements  as  scheduled.  Denote  by  4/ (4+1  )  the  desired  error 
covariance  of  the  state  estimation  at  time  4+1 .  Then  the  need  of  a  sensor  measure¬ 
ment  at  4+i  for  this  target,  according  to  the  covariance  control  method,  is 

«(4+i)  =  -min(eig{/Jd(4+1)  -  Pk+i\k})  (12.34) 

where  the  negative  sign  has  the  following  implication:  A  positive  value  of  the  eigen¬ 
value  difference  implies  that  the  desired  covariance  requirement  is  not  met.  The  goal 
of  covariance  control  is  to  minimize  the  total  sensing  cost  so  that  all  the  desired  co- 
variance  requirements  are  met.  If  we  use  the  same  notation  as  the  cost  for  the 
observer  i  to  sense  the  target  j  at  4+1,  then  the  covariance  control  tries  to  solve  the 
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following  optimization  problem. 


min  cijXij  (12.35) 

Xij 

subject  to  nij(tk+ 1)  =  -  mm(eig{ Pdij(tk+i)  -  /Vnifc+lCXi/)})  <  0, 

i  =  l, ...,  M,  j  =  \,  ...,N  (12.36) 

and  Xij  e  {0,  1}.  In  this  setting,  more  than  one  observer  can  sense  the  same  target 
at  the  same  time.  The  optimization  problem  is  combinatorial  in  nature  and  one  has 
to  evaluate  2MN  possible  sensor-to-target  combinations  in  general,  which  is  compu¬ 
tationally  prohibitive.  Alternatively,  a  suboptimal  need-based  greedy  algorithm  has 
been  proposed  [10],  It  also  considers  the  case  where  certain  constraints  can  not  be 
met,  i.e.,  the  desired  covariance  is  unachievable  even  with  all  the  sensing  resources. 


12.4.3  Game  Theoretic  Covariance  Prediction  for  Sensor 
Management 

In  the  formulation  of  the  sensor  management  problem,  we  need  the  filter  to  provide 
the  information  on  the  state  estimation  error  covariance  which  is  based  on  the  or¬ 
bital  trajectory  propagation  assuming  that  the  target  does  not  maneuver.  If  the  target 
maneuvers,  then  the  filter  calculated  error  covariance  assuming  the  non-maneuver 
motion  model  will  be  too  optimistic.  There  are  two  possible  approaches  to  account 
for  the  target  maneuver  motion.  One  is  to  detect  target  maneuver  and  estimate  its  on¬ 
set  time  as  quickly  as  possible  [15].  Then  the  filter  will  be  adjusted  with  larger  pro¬ 
cess  noise  covariance  to  account  for  the  target  maneuvering  motion.  Alternatively, 
one  can  design  a  few  typical  target  maneuvering  motion  models  and  run  a  mul¬ 
tiple  model  estimator  with  both  non- maneuver  and  maneuver  motion  models  [2], 
The  multiple  model  filter  will  provide  the  model  conditioned  state  estimation  error 
covariances  as  well  as  the  unconditional  error  covariance  for  sensor  management 
purposes.  Note  that  the  multiple  model  estimator  does  not  make  a  hard  decision 
on  which  target  motion  model  is  in  effect  at  any  particular  time,  but  evaluates  the 
probability  of  each  model.  The  corresponding  unconditional  covariance  immedi¬ 
ately  after  target  maneuver  onset  time  can  still  be  very  optimistic,  which  is  needed 
to  support  the  evidence  that  a  maneuvering  motion  model  is  more  likely  than  a  non¬ 
maneuvering  one.  As  a  consequence,  the  scheduled  sensing  action  in  response  to 
the  target  maneuver  based  on  the  unconditional  covariance  from  a  multiple  model 
estimator  can  be  too  late  for  evasive  target  motion. 

We  propose  to  use  generalized  Page’s  test  for  detecting  target  maneuver  [15]  and 
apply  the  model  conditioned  error  covariance  from  each  filter  in  the  sensor  manage¬ 
ment.  Denoted  by  Sm{tk+\)  the  set  of  targets  being  classified  as  in  the  maneuvering 
mode  and  S-m(tk+ 1)  the  set  of  targets  in  the  nonmaneuvering  model,  respectively. 
We  apply  covariance  control  for  sensor-to-target  allocation  only  to  those  targets  in 
Sm (tk+ 1 )  and  use  the  remaining  sensing  resources  to  those  targets  in  S-m(tk+ 1)  by 
maximizing  the  information  gain.  The  optimization  problem  becomes 
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min  CtjXij 

Xij 

(12.37) 

subject  to 

«//(4+i)  =  -min(eig {Pdij(tk+\)  ~  Pk+i\k+i(Xij)})  <  0, 

/  =  1 ,  . . . ,  M ,  j  G  Sm  1 ) 

(12.38) 

M 

1  ’  J  ^  S—m(tk+ 1) 

i=  1 

(12.39) 

§ 

ii 

VI 

>< 

Co 

u; 

(12.40) 

and  Xij  £  {0,  1}.  The  cost  Cjj  is 


c  _  |  Tr{//(xy((ffc+i)))  Rt  (tk+i)Hi  (xj  ((4+1 ))) }  j  C  S—m  (4+ 1 )  Q2  41) 
lJ  jo  j  e  Sm(tk+\) 

Given  that  target  j  e  S-m(tk),  we  assume  that  observer  i  will  declare  j  C  Sm  (tk+ 1 ) 
if  the  predicted  target  location  has  Earth  blockage  to  the  observer  i  at  tk+ \.  If  the 
target  is  declared  as  in  the  maneuvering  mode,  then  the  predicted  error  covariance 
will  be  based  on  the  worst  case  scenario  of  the  pursuit-evasion  game  between  the 
observer  and  the  target  (see  Appendix  4  for  details). 


12.5  Simulation  Study 
12.5.1  Scenario  Description 

We  consider  a  small  scale  space  target  tracking  scenario  where  four  satellite  ob¬ 
servers  collaboratively  track  two  satellite  targets.  The  nominal  orbital  trajectories 
are  generated  from  realistic  satellite  targets  selected  in  the  SpaceTrack  database. 
The  four  observer  satellites  are  (1)  ARIANE  44L,  (2)  OPS  856,  (3)  VANGUARD  1, 
and  (4)  ECHO  1.  They  are  in  low-Earth  orbits.  The  two  target  satellites  are: 
(1)  ECHOSTAR  10,  and  (2)  COSMOS  2350.  They  are  in  geostationary  orbits.  The 
simulation  was  based  on  the  software  package  that  utilizes  the  general  mission  anal¬ 
ysis  tool  (GMAT).1  Figure  12.1  shows  the  orbital  trajectories  of  the  observers  and 
targets.  The  associated  tracking  errors  were  obtained  based  on  the  recursive  linear 
minimum  mean  square  error  filter  when  sensors  are  assigned  to  targets  according  to 
the  nonmaneuvering  motion.  The  error  increases  in  some  time  segments  are  due  to 
the  Earth  blockage  where  no  measurements  are  available  from  any  observer. 


'General  Mission  Analysis  Tool,  released  by  National  Aeronautics  and  Space  Administration 
(NASA),  available  at  http://gmat.gsfc.nasa.gov/. 
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Fig.  12.1  The  space  target  tracking  scenario  where  four  observers  collaboratively  track  two  targets 

12.5.2  Performance  Comparison 

We  now  consider  the  case  in  which  the  target  performs  an  unknown  thrust  maneu¬ 
ver  that  changes  the  eccentricity  of  its  orbit.  In  particular,  both  targets  are  initially 
in  the  GEO  orbit  and  at  time  t  —  1000  s,  target  1  performs  a  1  s  burn  that  produces  a 
specific  thrust,  i.e.,  an  acceleration  w  =  (0  0.3  0]r  km/s2,  while,  at  time  t  —  1500  s, 
target  2  performs  a  1  s  burn  that  produces  a  specific  thrust  w  =  [0  0.5  0] 7  km/s2. 
The  eccentricity  change  of  target  1  after  the  burn  is  around  e  0.35  while  the  eccen¬ 
tricity  change  of  target  2  after  the  second  bum  is  e  0.59.  Another  type  of  target 
maneuver  is  the  inclination  change  produced  by  a  specific  thrust.  Two  inclination 
changes  are  generated  with  i  ~  0.16  for  target  1  and  i  ~  0.09  for  target  2.  Each 
observer  has  the  minimal  sampling  interval  of  50  s  and  we  assume  that  all  observers 
are  synchronized  and  the  sensor  manager  can  make  centralized  coordination  based 
on  the  centralized  estimator  for  each  target.  We  consider  two  cases:  (i)  the  observer 
has  range  and  angle  measurements  with  standard  deviations  0.1  km  and  10  mrad,  re¬ 
spectively,  and  (ii)  the  observer  has  range,  angle  and  range  rate  measurements  with 
standard  deviations  0.1  km,  10  mrad,  2  m/s,  respectively. 

We  applied  the  generalized  Page’s  test  (GPT)  without  and  with  range  rate  mea¬ 
surement  for  maneuver  detection  while  the  filter  update  of  state  estimate  does  not 
use  the  range  rate  measurement  [15].  The  reason  is  that  the  nonlinear  filter  designed 
assuming  non-maneuver  target  motion  is  sensitive  to  the  model  mismatch  in  the 
range  rate  when  target  maneuvers.  The  thresholds  of  the  generalized  Page’s  test  for 
both  case  (i)  and  case  (ii)  were  chosen  to  have  the  false  alarm  probability  Pfa  =  1%. 
The  average  delays  of  target  maneuver  onset  detection  for  both  cases  are  measured 
in  terms  of  the  average  number  of  observations  from  the  maneuver  onset  time  to 
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Table  12.1  Comparison  of  tracking  accuracy  for  various  orbital  changes  of  the  targets 


Cases 

e  -»■  0.35 

e  -*  0.59 

i  -*  0.09 

(  -*  0.16 

(i)  Average  delay  per  observations 

7.3 

6.4 

12.1 

9.6 

(i)  PTR,  peak  position  error  (km) 

33.7 

54.8 

14.5 

18.9 

(i)  IMM,  peak  position  error  (km) 

48.1 

66.5 

12.4 

22.3 

(i)  PTR,  peak  velocity  error  (km/s) 

0.38 

0.40 

0.22 

0.25 

(i)  IMM,  peak  velocity  error  (km/s) 

0.41 

0.44 

0.21 

0.26 

(ii )  Average  delay  per  observations 

1.3 

1.0 

2.9 

2.4 

(ii )  PTR,  peak  position  error  (km) 

6.9 

7.4 

3.1 

3.8 

(ii)  IMM,  peak  position  error  (km) 

8.3 

11.2 

3.2 

4.1 

(ii )  PTR,  peak  velocity  error  (km/s ) 

0.28 

0.31 

0.16 

0.19 

(ii)  IMM,  peak  velocity  error  (km/s) 

0.30 

0.34 

0.16 

0.21 

the  declaration  of  the  target  maneuver.  Once  target  maneuver  is  declared,  then  an¬ 
other  filter  assuming  the  white  noise  acceleration  with  process  noise  spectrum  of 
0.6  km/s2  was  used  along  both  tangential  and  normal  directions  of  the  estimated  tar¬ 
get  motion.  Alternatively,  an  interacting  multiple  model  (1MM)  estimator  [2]  with 
nonmaneuver  and  maneuver  motion  models  using  the  same  parameter  settings  as 
the  model  switching  filters  embedded  in  the  target  maneuvering  detector  was  used 
to  compare  the  tracking  accuracy.  Table  12.1  compares  the  peak  errors  in  position 
and  velocity  for  each  target  maneuvering  motion  using  model  switching  filter  and 
the  IMM  estimator.  The  average  detection  delays  for  both  cases  are  also  shown 
in  Table  12.1  for  GPT  algorithm.  We  can  see  that  the  range  rate  measurement,  if 
available,  can  improve  the  average  detection  delay  of  target  maneuver  significantly 
and  thus  reducing  the  peak  estimation  errors  in  both  target  position  and  velocity. 
The  model  switching  filter  using  GPT  has  better  tracking  accuracy  than  the  IMM 
estimator  in  most  cases  even  for  the  peak  errors.  It  should  be  clear  that  the  filter 
based  on  nonmaneuver  motion  model  outperforms  the  IMM  estimator  during  the 
segment  that  the  target  does  not  have  an  orbital  change.  Interestingly,  even  though 
the  inclination  change  takes  longer  time  to  detect  compared  with  the  eccentricity 
change  for  both  targets,  the  resulting  peak  estimation  errors  in  position  and  velocity 
are  relatively  smaller  for  both  the  model  switching  filter  and  the  IMM  estimator. 
The  extensive  comparison  among  other  nonlinear  filtering  methods  for  space  target 
tracking  can  be  found  in  [4]. 

Next,  we  assume  that  both  target  1  and  target  2  will  choose  their  maneuver  on¬ 
set  times  intelligently  based  on  their  geometries  to  the  observers.  Both  targets  can 
have  a  maximum  acceleration  of  0.05  km/s2  with  a  maximum  of  10  s  burn.  We  as¬ 
sume  that  each  observer  can  measure  target  range,  angle,  and  range  rate  with  the 
same  accuracy  as  in  the  case  (ii)  of  the  tracking  scenario  considered  previously.  We 
implemented  the  following  configurations  of  the  sensor  management  schemes  to 
determine  which  observer  measures  which  target  at  a  certain  time  instance,  (i)  In¬ 
formation  based  method:  Sensors  are  allocated  with  a  uniform  sampling  interval  of 
50  s  to  maximize  the  information  gain,  (ii)  Covariance  control  based  method:  Sen- 
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Table  12.2  Performance  comparison  for  various  sensor  management  methods 

Configuration  Peak  position  Average  position  Peak  velocity  Average  velocity  Average  sampling 
error  (km)  error  (km)  error  (krn/s)  error  (km/s)  interval  (s) 


(i)  Target  1 

89.7 

24.4 

2.1 

0.16 

50 

(i)  Target  2 

76.3 

14.2 

1.6 

0.14 

50 

(ii)  Target  1 

16.7 

2.4 

0.36 

0.10 

13.5 

(ii)  Target  2 

15.2 

1.8 

0.29 

0.09 

16.8 

(iii)  Target  1 

12.2 

1.8 

0.28 

0.08 

37.6 

(iii)  Target  2 

10.9 

1.3 

0.27 

0.08 

39.2 

sors  are  scheduled  with  the  sampling  interval  such  that  the  desired  position  error  is 
within  10  km  for  each  target.  The  state  prediction  error  covariance  for  nonmaneuver 
or  maneuver  motion  model  is  computed  based  on  the  posterior  Cramer-Rao  lower 
bound,  (iii)  Game  theoretic  method:  Sensors  are  allocated  by  maximizing  the  infor¬ 
mation  gain  for  nonmaneuvering  targets  and  covariance  control  will  be  applied  to 
maneuvering  targets.  The  maneuvering  onset  time  is  predicted  based  on  the  pursuit 
evasion  game  for  each  observer-to-target  pairing.  In  all  three  configurations,  model 
switching  filter  is  used  for  tracking  each  target.  We  compare  the  peak  and  average 
errors  in  position  and  velocity  for  both  targets  as  well  as  the  average  sampling  inter¬ 
val  from  all  observers  for  configurations  (i)-(iii).  The  results  are  listed  in  Table  12.2. 
We  can  see  that  the  information  based  method  (configuration  (i))  yields  the  largest 
peak  errors  among  the  three  schemes.  This  is  due  to  the  fact  that  both  targets  apply 
evasive  maneuver  so  that  the  peak  errors  will  be  much  larger  compared  with  the  re¬ 
sults  in  Table  12. 1  based  on  the  same  tracker  design  and  sampling  interval.  In  order 
to  achieve  the  desired  position  error,  covariance  control  method  (configuration  (ii)) 
achieves  much  smaller  peak  and  average  errors  compared  with  configuration  (i)  at 
the  price  of  making  the  sampling  interval  much  shorter.  Note  that  the  peak  position 
errors  are  larger  than  10  km  for  both  targets  owing  to  the  detection  delay  of  maneu¬ 
vering  onset  time.  The  proposed  game  theoretic  method  (configuration  (iii))  has  the 
smallest  peak  and  average  errors  because  of  the  prediction  of  maneuvering  onset 
time  by  modeling  target  evasive  motion  from  the  pursuit  evasion  game.  Note  that 
it  also  has  longer  average  sampling  interval  than  that  of  configuration  (ii)  due  to  a 
quicker  transient  when  each  target  stops  its  burn,  indicating  possible  energy  saving 
for  the  overall  system.  The  sensing  resources  saved  with  configuration  (iii)  can  also 
be  applied  to  search  other  potential  targets  in  some  designated  cells.  One  possible 
approach  to  perform  joint  search  and  tracking  of  space  targets  was  discussed  in  [5], 


12.6  Summary  and  Conclusions 

In  space  target  tracking  by  satellite  observers,  it  is  crucial  to  assign  the  appropriate 
sensor  set  to  each  target  and  minimizes  the  Earth  blockage  period.  With  complete 
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knowledge  of  the  space  borne  observers,  a  target  may  engage  its  evasive  maneu¬ 
vering  motion  immediately  after  the  Earth  blockage  occurs  and  change  its  orbit 
to  maximize  the  duration  of  the  Earth  blockage.  We  presented  a  game  theoretic 
model  for  the  determination  of  maneuvering  onset  time  and  consequently,  the  co- 
variance  control  is  applied  to  the  maneuvering  targets  in  the  sensor-to-target  assign¬ 
ment.  For  nonmaneuvering  targets,  we  try  to  maximize  the  total  information  gain 
by  selecting  sensors  that  will  provide  the  most  informative  measurements  on  the 
target’s  state.  We  simulated  a  multi-observer  multi-target  tracking  scenario  where 
four  LEO  observers  collaboratively  track  two  GEO  targets.  We  found  that  the  mul¬ 
tiple  model  estimator  assuming  random  maneuvering  onset  time  yields  much  larger 
estimation  error  compared  with  the  model  switching  tracker  based  on  maneuver¬ 
ing  detection  from  the  solution  to  the  pursuit-evasion  game.  In  addition,  sensor 
assignment  based  on  maximum  information  gain  can  lead  to  large  tracking  error 
for  evasive  targets  while  using  the  same  desired  error  covariance  for  all  targets  can 
only  alleviate  the  issue  at  the  price  of  more  frequent  revisit  time  for  each  target. 
Fortunately,  the  sensor  assignment  based  on  covariance  control  for  maneuvering 
targets  and  maximum  information  gain  for  nonmaneuvering  targets  achieves  a  rea¬ 
sonable  tradeoff  between  the  tracking  accuracy  and  the  consumption  of  sensing  re¬ 
sources. 

There  are  many  avenues  to  extend  the  existing  work  in  order  to  achieve  space 
situational  awareness.  First,  target  intent  can  be  inferred  based  on  its  orbital  history. 
It  is  of  great  interest  to  separate  the  nonevasive  and  evasive  orbital  maneuvering 
motions  and  allocate  the  sensing  resources  accordingly.  Second,  our  model  of  the 
pursuit-evasion  game  relies  on  the  complete  knowledge  of  the  observer’s  and  tar¬ 
get’s  state  which  may  not  be  known  to  both  players  in  real  life.  This  poses  challenges 
in  the  development  of  the  game  theoretic  model  with  incomplete  information  which 
is  computationally  tractable  for  the  sensor  management  to  allocate  sensing  resources 
ahead  of  time.  Finally,  the  current  nonlinear  filter  does  not  consider  the  clutter  and 
closely  spaced  targets  where  imperfect  data  association  has  to  be  handled  by  the  fil¬ 
tering  algorithm.  It  should  be  noted  that  the  posterior  Cramer-Rao  lower  bound  for 
single  target  tracking  with  random  clutter  and  imperfect  detection  [19]  can  be  read¬ 
ily  applied  to  the  covariance  control.  However,  it  is  still  an  open  research  problem 
to  design  efficient  nonlinear  filtering  method  that  can  achieve  the  theoretical  bound 
of  the  estimation  error  covariance. 
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Appendix  1:  Conversion  of  the  Coordinate  Systems 

The  following  conversion  schemes  among  different  coordinate  systems  are  based 
on  [3],  Given  the  position  r  =  [f  r]  £]'  in  the  ECEF  frame,  the  latitude  (p,  longitude 
X  and  altitude  h  (which  are  the  three  components  of  rgeo)  are  determined  by 


1 2  Orbital  Evasive  Target  Tracking  and  Sensor  Management 


249 


STEP  1  <p  =  0 
repeat 

tPo\d  —  V 

Dy  =  Re[\  —  €e  sin2  <p0id]-' 

't,  +  D^e  singpid 

rj2 


STEP  2 


w  =  atan .  , 

until  |  <p  —  ^oidl  <TOL 
X  —  atan(|) 


h  — 


c 


stn^> 


M  l~€2e) 


where  Re  =  6378.137  km  and  ee  =  0.0818191  are  the  equatorial  radius  and  eccen¬ 
tricity  of  the  Earth,  respectively.  TOL  is  the  error  tolerance  (e.g.,  10-10)  and  the 
convergence  occurs  normally  within  10  iterations. 

The  origin  of  the  local  Cartesian  frame  O  is  given  by 


0=  [0  Dve^  sincpcoscp  sin2<p  —  l)]/ 

and  the  rotation  matrix  is  given  by 


.4  = 


—  sink 
—  simpcosk 
cos  cp  cos  X 


cosk  0 

—  sin  <p  sin  k  cos  q> 
cos  (p  sink  sin</> 


(12.42) 


(12.43) 


The  position  r  in  the  local  Cartesian  frame  is  given  by 


rioc  =  Ar  +  O 


(12.44) 


Appendix  2:  Keplerian  Elements 


The  specific  angular  momentum  lies  normal  to  the  orbital  plane  given  by  h  =  r  x  v 
with  magnitude  h  =  ||h||.  Inclination  is  the  angle  between  the  equatorial  plane  and 

the  orbital  plane  given  by  i  =  cos- 1  (jf  )  where  h-  is  the  "-component  of  h.  Eccen¬ 
tricity  of  the  orbit  is  given  by 


A 

e  = 


(12.45) 


with  magnitude  e  —  ||e||.  The  longitude  of  the  ascending  node  is  given  by 


£2  = 


cos  !(^f) 
2;r-cos-1(^) 


ny  >  0 
ny  <  0 


(12.46) 
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where  n  is  the  vector  pointing  towards  the  ascending  node  with  magnitude  n  —  ||n|| . 
The  argument  of  perigee  is  angle  between  the  node  line  and  the  eccentricity  vector 
given  by 


A 

O)  = 


cos  1  (^f )  ez>0 

27T-cos_1(^)  ez  <  0 


with  the  convention  that 


co  —  cos 


(12.47) 


(12.48) 


for  an  equatorial  orbit.  The  true  anomaly  v  is  the  angle  between  the  eccentricity 
vector  and  the  target’s  position  vector  given  by 


A  |  COS  J(P)  Vr  >  0 

I  27T  -  COS_1(^)  Vr  <0 


(12.49) 


with  the  convention  that  v  =  cos  1  (y)  for  a  circular  orbit.  The  eccentric  anomaly 
is  the  angle  between  apogee  and  the  current  position  of  the  target  given  by 


E  —  cos 


1  —  r  /a 


(12.50) 


where  a  is  the  orbit’s  semi-major  axis  given  by  a  —  1  2  .  The  mean  anomaly  is 


M  =  E  —  e  sin  E .  The  orbital  period  is  given  by  T  =  2tt  J  . 

The  six  Keplerian  elements  are  [a,  i.  Q ,  co,  e.  M\.  The  orbit  of  a  space  target  can 
be  fully  determined  by  the  parameter  set  {j,  Q ,  co,  T .  e,  M }  with  initial  condition 
given  by  the  target  position  at  any  particular  time  [17].  The  angles  {i,  Q .  co}  trans¬ 
form  the  inertial  frame  to  the  orbital  frame  while  T  and  e  specify  the  size  and  shape 
of  the  ellipsoidal  orbit.  The  time  dependent  parameter  v(t)  represents  the  position 
of  the  target  along  its  orbit  in  the  polar  coordinate  system. 


Appendix  3:  Algorithm  for  Orbital  State  Propagation 

Presented  below  is  an  algorithm  that  propagates  the  state  of  an  object  in  an  or¬ 
bital  trajectory  around  the  Earth  following  [3],  Both  the  trajectory  propagation 
and  the  corresponding  Jacobian  matrix  of  the  nonlinear  orbital  equation  are  given. 
Let  x{t)'  —  [r(f)'  f(f)']  be  the  unknown  state  to  be  computed  at  time  t ,  given 
the  state  =  x(to)/  =  [fo^o]  at  the  time  tp.  The  gravitational  parameter  /x  = 
3.986012  x  105  km3 /sec2  and  the  convergence  check  parameter  TOL  =  10-10  are 
used. 
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STEP  1 


STEP  2 

STEP  3 

STEP  4 


STEP  5 


ro  :=  II r0||;  vo:=l|rol|;  qo  '■=  -r0r0 

P 

2  115  1  -  a0r0 

«o  := - ;  PO  ■=  — 

ro  P  V/1 

ao(t-to)  0  2 

o'  := - ;  p  := 

V/1 

1  -  cos(V?) .  V?  _  sin(ViS) 

0  ;  5:__ 

3  i  2  ,  rQ 

x  pool  s  +  qoa  c  H - a 

Vp 


ui  2  A0 

—  :=  poa~c  +  qo<x(l  -  sfi)  +  —— 
dot 

[(t-to)-x] 
if([(f-f0)-r]  >TOL) 


a  :=  a  + 


dr 

do? 


STEP  6 


STEP  7 


gotoSTEP  3 


2 

O'  c 


3 

as 


/  :=  1 - ;  g '■=(!- to) - — 

ro  y/P 

r(t)  :=  fr0  +  giro',  r:=||r(f)|| 


\  ro 


fo8-l)  -  ;  g:=  1 

r 


OL  C 

- ;  r(f)  := /r0  +  gr0 


The  above  steps  yield  the  required  state  x(t)'  =  [rtf/rtf )']  at  time  t .  In  order  to 
predict  the  covariance  of  the  position  rtf )  by  propagating  the  covariance  of  r(fo) 
from  to  to  f,  we  need  to  compute  the  6x3  matrix  Vlor(f ).  The  computation  of  this 
matrix  involves  the  following  additional  steps. 


STEP  8 


V.,0r0  := 


ro 

0 


Vx0qo 


V.^flo  :=  (V.vor0) 


+  (Y^ro) 


VxoPo  :=  (VT0ro) 


/  ~ao\ 

\Vp) 


(V.voflo) 


(5) 
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STEP  9 

STEP  10 


ds  c  —  3s  dc  1  —  s/3  —  2c 
d/3  _  2/S  d/T~  2/3 


b\  :=  (Yv0^o) (—a2c)  +  (VAo/?0)(-Q!3.j)  +  (VXor0) 

bi  :=  (YY0fl0)(-a2) 


STEP  11 


A:= 


3p0(xzs +  2q0ac+ ^  Pool3  ^  +  q0az  ^ 
2«o«  — 1 


„2  dr 


[(VY0a)(VYoJ6)]  :=  [bxb2]A 


-l 


V,0  /  := 


(VA-0r0) 


V  ro 


(Vxo«)(2c)-(VYoyS)(a0 


/  dj\l 

T—  a2~\ 

v-v0g  := 

(VY0«)(3s)-(VY0 /))(«— J 

.  V/1, 

step  12  VAQr(r) 


//3 


+  (VA0/)r'+(VA-0g)r'. 


Appendix  4:  Pursuit  Evasion  Game  in  a  2D  Plane 


Consider  the  space  target  orbiting  the  Earth  with  the  polar  coordinate  system  fixed 
on  the  Earth’s  center.  The  motion  of  the  target  is  given  by 


r  -  rd2 
rd  +  2  r6 


F  sin  a 
m 


F  cos  a 


m 


(12.51) 

(12.52) 


where  a  is  the  angle  the  thrust  vector  and  the  local  horizontal  as  shown  in  Fig.  12.2. 
Denote  by  vq  and  vr  the  tangential  and  radial  velocities  of  the  target,  respectively. 
The  dynamic  equation  of  the  space  target  can  be  written  as 


i)r  - 


r 


.  ,  VrVg 

ve  H - 

r 


F  sin  a 


m 


F  cos  a 


m 


(12.53) 

(12.54) 


It  is  desirable  to  normalize  the  parameters  with  respect  to  a  reference  circular 
orbit  with  radius  ro  and  velocity  no,  so  that  significant  figures  will  not  be  lost  due  to 
linearizing  the  target  motion  equation.  Define  the  following  normalized  state  vari¬ 
ables. 

r 

x\  =  — 

ro 


(12.55) 
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Fig.  12.2  Simplified  2D 
geometry  of  a  space  target 


Vr 

*2  —  ~ 
VO 

ve 

*3  =  — 

vo 

X4  =  6 


(12.56) 

(12.57) 

(12.58) 


Let  r  =  y-t.  The  state  equation  with  respect  to  r  can  be  written  as 

x\  —  X2  (12.59) 

x 2  i 

xi  —  — - y  +  Tosina’  (12.60) 

xi  xi 

where  Fq  —  is  a  constant  depending  the  thrust  F  and  target’s  mass  in. 


X1X3 

X3  = - b  Fq  cos  a 

x\ 

x3 

X4  =  — 

XI 


(12.61) 

(12.62) 


In  the  standard  pursuit  evasion  game,  the  objective  is  to  find  the  minimax  solu¬ 
tion,  if  exists,  to  the  objective  function 


J  =  <P{x(tf))  (12.63) 

with  the  state  dynamics  given  by 

x=/(x,  u,v,t)  (12.64) 

and  the  initial  condition  x(fo)  =  xo  as  well  as  the  terminal  constraint  1 {r(x(tf))  =  0. 
Here  u  and  v  represent  the  controls  associated  with  the  pursuer  and  the  evader, 
respectively.  The  goal  is  to  determine  {m*,  i>*}  such  that 

j(u* ,  v)  <  J(u* ,  r>*)  <  J(u,  n*) 


(12.65) 
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The  necessary  condition  for  the  minimax  solution  to  exist  is  that  the  costate  A  satis¬ 
fies  the  following  transversality  condition. 

k  =  +  (12.66) 

H(tf)  =  4>,(tf)  +  v^t(tf)  (12.67) 

where  v  is  a  Lagrange  multiplier  and  H  is  the  Hamiltonian  associated  with  A  that 
has  to  be  optimized 

H*  =  maxminH(x,  A,  u,  v,  t)  (12.68) 

V  U 

It  has  been  shown  in  [9]  that  at  the  optimal  solution,  the  thrust  angles  of  the  pursuer 
and  the  evader  are  the  same. 
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